# connectivity.R
# building "connectivity" matrix C with periodic boundary conditions

periodic.matrix <- function(n)
{
	stopifnot(n >= 3)

	A <- diag(-8,n)
	A[row(A)==col(A)+1] <- rep(1,n-1)
	A[row(A)==col(A)-1] <- rep(1,n-1)
	A[1,n] <- 1
	A[n,1] <- 1

	B <- A
	diag(B) <- 1

	C <- matrix(0,n^2,n^2)
	for (i in 1:n) {
		C[(1+n*(i-1)):(n*(i)),(1+n*(i-1)):(n*(i))] <- A
		if (i!=n) {
			C[(1+n*i):(n*(i+1)),(1+n*(i-1)):(n*(i))] <- B
			C[(1+n*(i-1)):(n*(i)),(1+n*i):(n*(i+1))] <- B
		}
	}
	C[1:n,(1+n*(n-1)):(n*n)] <- B
	C[(1+n*(n-1)):(n*n),1:n] <- B
	C
}
